Counterflow superfluid of polaron pairs in Bose-Fermi mixtures in optical lattices 
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We study the quantum phases of one-dimensional Bose-Fermi mixtures in optical lattices. Assum- 
ing repulsive interparticle interactions, equal mass, and unit total filling, we calculate the ground- 
state phase diagram by means of both Tomonaga-Luttinger liquid theory and time-evolving block 
decimation method. We demonstrate the existence of a counterflow superfluid (CFSF) phase of 
polaron pairs, which are composite particles consisting of two fermions and two bosonic holes, in 
a broad range of the parameter space. We find that this phase naturally emerges in 174 Yb- 173 Yb 
mixtures, realized in recent experiments, at low temperatures. 

PACS numbers: 03.75.-b, 67.85.Pq, 05.30.Rt 
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The unprecedented control that has been achieved in 
ultra-cold atom systems in optical lattices has gener- 
ated a new frontier in exploring quantum phases in the 
strongly correlated regime. Novel quantum phases have 
been created in mixtures of different hyperfine states 0- 
13] , atomic species 0-0] , and isotopes d, Q in optical lat- 
tices. In these experiments, a wide range of features, such 
as the statistics of particles, the mass and density ratios, 
and interparticle interactions, can be precisely varied. 
This, in turn, has led to the prediction of numerous fur- 
ther exotic phases that can be studied in these systems, 
UHj], paired sup erfluids fl3l — 



including supersolids 
HH, counterflow superfluids (CFSF) 



Tomonaga-Luttinger liquids (TLL) of polarons 




and 



Recently, experiments on a Bose-Fermi (BF) mixture 
of Ytterbium isotopes ( 174 Yb- 173 Yb) in optical lattices 
and the realization of a novel dual Mott insulator have 
been reported 0. In this phase the total density (n^j + 
Tif,i), where nuf),i is the bosonic (fermionic) density at 
lattice site i, is pinned to unit filling, while the individual 
densities (nf,,t) and (n/,j) arc unrestricted and assumed 
to be close to half-filling. This suggests that the degree 



of freedom that has to be considered is ct ,c& , 



where i 



(cb,i) is the fermionic creation (bosonic annihilation) op- 
erator at site i, because the creation of one particle type is 
matched with the annihilation of the other, thus keeping 
the density at unity. Thus, the sy_stem can be regarded as 
a liquid of composite fermions [13, EH , which can also be 



considered as a limiting case of polarons ll| . If this mix- 
ture was either a Bosc-Bose (BB) or Fermi-Fermi (FF) 
mixture, these composite particles are bosonic. For a BB 
mixture, say, these bosonic particle-hole pairs can con- 
dense into a CFSF phase [l4, 16- 18]. However, the na- 



ture of BF mixtures is fundamentally distinct, because 
the particle-hole pairs are fermionic, and thus cannot 
function as an order parameter. Rather, these compos- 
ite fermions have to form pairs themselves to condense. 
This constitutes a Bose-Einstein condensation of quar- 
tets consisting of two fermions and two bosonic holes, 
thus creating a polaron-paired CFSF (PP-CFSF) phase. 



Previous studies of these systems in the strong-coupling 
limit have predicted phase separation (PS), Fermi liquid 
of polarons, spin-density wave (SDW), and CFSF with p- 
wave pairing of the polarons \yA [l8( ] . This p-wave CFSF 
is fascinating also in the sense that it can be regarded as 
a topological superfluid , which has been extensively 
discussed in the contexts of liquid 3 He [Hf and the com- 
pound of Sr2RuC>4 [2^]. However, in direct simulations of 
the BF- Hubbard model [H [27rl3^. the presence of this 
exotic state has not been confirmed. 

In this Letter, we show that the PP-CFSF phase is 
naturally realized in Yb mixtures, which constitutes the 
first numerical demonstration of this phase in the BF- 
Hubbard model. In Fig.Q]wc show the quantum phase di- 
agram of the one-dimensional (ID) BF-Hubbard model, 
obtained with the quasi-exact time-evolving block deci- 
mation (TEBD) method (HI, supported by TLL theory. 
The Hamiltonian of the system is [34| : 

H = - *«(4jc»j+i + H - c -) + e J n ».j 



« 3 



a 3 



where t a is the tunneling energy of particle type a = b, f. 



€j denotes the external potential, and n a , 



- J 



and Ubb denote the on-site interactions. For a homoge- 
neous system, we define the filling fractions v a = {n a ,j}- 
Inspired by 174 Yb- 173 Yb mixtures reported in Ref. 0, we 
assume that U aa > > 0, t = tf, = tf, and v = vt + Vf, = 1. 
In the phase diagrams in Fig. [TJ we vary t/Ubb, Ubf /Ubb, 
and Vb — Vf. We find that the PP-CFSF phase occupies 
a broad regime between the SDW and PS regions, which 
is approximately 1 < Ubb /Ubf < 2 and t/Ubb < 0.3 when 
v b = v f = 0.5. Since U bf /U bb ~ 1.32 in 174 Yb- 173 Yb 
mixtures, the ground state of the system is expected to 
be the PP-CFSF state. In the following, we first map 
out the ground-state phase diagram of the homogeneous 
system (ej — 0) varying t = t/Ubb, u = Ubf/Ubb, and 
i\v = v b - Vf to find the PP-CFSF phase. To con- 
nect closely to experiments, we also confirm that the PP- 
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FIG. 1: (color online) Ground-state phase diagrams for the Bose-Fermi Hubbard model of Eq. |T]) with equal hoppings, i.e., 
t = tb = tf. The emerging phases include bosonic superfluid (b-SF), spin-density wave (SDW), polaron-paired counterflow 
superfluid (PP-CFSF), and phase separation (PS). The phase diagrams from TLL theory (a) and TEBD (b) are depicted 
in the (Ubf /Ubb, t/t/i,i>)-plane, where Vb = Vf — 0.5. The purple thick-solid line represents the phase boundary between the 
two-component TLL and the Mott insulator, and the strong-coupling regime is achieved below the dashed line. The purple thin- 
solid lines separate the phases in the Mott insulating region. The green dot-dashed line in (a) separates the fermionic QLROs, 
i.e. polaron pairing (PP) and density wave (DW), in the two-component TLL region. The black dotted line in (b) represents 
Ubf /Ubb = 1.32, corresponding to the 174 Yb- 173 Yb mixture of Ref. 0. Our phase boundary is consistent with the results of 
Ref. [29^ | . shown as an orange dashed-dotted line, (c) The phase diagram from TEBD is depicted in the (Av, i/[4i,)-plane, with 
Ubf /Ubb = 1.4. We take L = 80 for the TEBD calculations. 



CFSF state emerges inside the Mott plateau with v = 1 
in the presence of a parabolic trapping potential. 

We first determine the phase diagram via TLL theory. 
We go to a continuum representation, Caj/d 1 ^ 2 — > ^) a (x), 
where d is the lattice constant, and express the particle 
operators through a bosonization identity [35l |36| : 



$ffh = [/>//*> + n //*>] 1/2 



E 

i odd/ even 



e imQ f/be i<$ f 



(2) 



where pf/ b = Vf/b/d are the average densities of the 
fcrmions/bosons, Hf/ b (x) are the low-fc parts (i.e. k <C 
l/pf,l/p b ) of the density fluctuations, $//{, arc the 
phase fields, and Qf/ b (x) = Kpf/ b x + 6f/ b (x), with 
9f/ b (x) = 7T J x dyllf / b (y). The right- and left-movers of 

the fermions are represented by 4>f,R/L ~ p 1 / 2 e ±lSf e l<Sf . 
The action describing the mixture is S = Sqj + So tb + 
Sbf + S u k- The terms Sa. a ar e 

1 f .0 / 1 



So.a 



d 2 r 



-(d T 6 a ) 2 + v a (d x a ) 2 ) (3) 



2wK a 

where v a and K a are the velocities and Luttinger param- 
eters of the fermionic and bosonic liquid, for U b f =0. We 
use H = 1 throughout the paper. The fermionic parame- 
ters are vt = 2tdsm(irvf) and Kf = 1, the bosonic ones 
are v b = Vf(l — 8ii/ b cos(irv b )) and K b = 1 + 8i sin(7rz^,)/7r 
assuming t <C 1. The interaction between bosons and 
fermions generates both linear terms, Sbf, and a non- 
linear term S u k, in the effective action. S b f is 

U bf f ;2„x a a a , V bf 



S b f 



d rd x 6\d x 02 



d z rd T 9 1 d T e 2 .(A) 



The second term is created during the renormalization 
group (RG) flow; its prefactor therefore has the initial 
value V b f(0) = 0. The non-linear Umklapp term S u k is 



S u k 



(27ra) 2 



d 2 r cos(26 1 + 26 2 ) 



(5) 



The renormalization group (RG) flow equations for this 
system were derived in Ref. its qualitative behav- 
ior is that of a Berezinskii-Kosterlitz-Thouless (BKT) 
transition (37j . S u k can cither be irrelevant, g u k 0, 
or relevant, g u k oo. In Fig. [TJa), the boundary 
between these regimes is depicted by the purple thick- 
solid curve. The relevant regime has two subrcgimes, 
namely the strong-coupling regime, in which g u k diverges 
rapidly, and the cross-over regime, in which g u k initially 
decreases. We estimate the boundary between these by 
the initial sign of the prefactor of the g u k flow equa- 
tion (see Ref. fllj). Using the above estimates we find 
U b f,c = 16t 2 /U bb (sin(irP b )) 2 (l - 4(tv b /U bb ) cos(nv b )), de- 
picted by the blue dashed line in Fig. QJa). The orders 
indicated in this figure are the dominant quasi-long-range 
order (QLRO), i.e. the order parameter 0(x) whose cor- 
relation function exhibits the slowest algebraic decay as 
(O T (x)O(0)) ~ \x\~ a , with a < 2, see e.g. Ref. [39j]. 
In the regime where S u k is irrelevant, we find that the 
most dominant QLRO is always the superfluid of bosons 
dressed with fermions whose order parameter is O db 
cxp(- 



-irj^> f)ip b with a real number r\ [22 
fermionic sector we find a competition of the 2kp compo- 
nent of the density operator of the fermions, fDW = pj, 
describing a spontaneous density modulation, and a po- 
laron pairing operator O = x/jf^ipf^exp^X^j,), with 
A a continuously varying parameter. 

For relevant S u k, the total density is frozen out, which 
suggests a strong-coupling expansion leading to the fol- 
lowing spinless Fermi-Hubbard model [13, EH , 



23J. For the 



£f = -jV(/t/ 3+ i+H.c.) + V 



j 



3 



(6) 



where U 



is the polaron annihilation oper- 



ator, consisting of an original fermion and a bosonic 
hole; m, = f]f h V = 2{t 2 + t 2 )/U bf - (U 2 )/U bb , and 



3 



(a) 





0.1 


» • 








1 


0.01 


• 


a 

o 








0.001 


■C B (r) 


1 




c 

u 


io-' 


•|C P (r)| 




• QpW 




I0" s 


t/ BF = 0.4 ! 



(b) 



•ICpMI . 

• C P p('') 

U Bf = 0.05 1 . 

1 2 5 TO 2(T 



FIG. 2: (color online) Blue squares, green diamonds, and red 
circles represent the correlation functions Cs(r), \Cp(r)\, and 
Cpp(r), respectively, for L = 80, Ubf /Ubb = 1.4, and £/J7m, = 
0.4 (a) and 0.05 (b). The plots are on a log- log scale. 



J = (2tf,t f) / Ub f ■ fj experiences no external field. The 
model of Eq. (|6|) is solvable by Bethe ansatz, see e.g. 
Ref. [iil: For V/J > 2 and v v = (rrij) = 0.5, i.e., 
(0 <)u < |(1 — g) 2 and Vf, = i>f, where g = tf/tb, 
Neel order with wave number k = ir/d emerges. This 
phase is not present when tb = tf. For V/J < —2, i.e. 
u > i(l + g) 2 , strong attraction leads to the collapse 
of the polaron gas, i.e. phase separation of the original 
mixture. The remaining region is a TLL of polarons. For 
positive V, i.e., u < ^(1 + g 2 ), the dominant QLRO is 
density wave ordering, corresponding to a SDW phase of 
the original mixture. For negative V, i.e. u > i(l + g 2 ), 
the dominant QLRO is triplet polaron-pairing, i.e. the 
PP-CFSF phase, as indicated in Fig.QJa). 

Having established this phase analytically, we now cor- 
roborate its existence numerically, and demonstrate its 
stability in a trapped system. We calculate the ground- 
state phase diagram using the TEBD method for open 
boundaries (33j via imaginary-time propagation. We fix 
the number of lattice sites to L = 80. 

The TEBD method is a variant of the density-matrix 
renormalization group (DMRG), which allows to accu- 
rately calculate ground states of much larger systems 
than those tractable with exact diagonalization used in 
Ref. [22j ■ In addition, it has the advantage that any cor- 
relation function can be calculated efficiently, in contrast 
to quantum Monte Carlo methods used in Refs. 28 - 30| . 
for example. Although DMRG has been applied to the 
BF-Hubbard model in Refs. (HSU, the PP-CFSF phase 



has not been identified and explored, as we do here. 

We first determine the transition from the two- 
component TLL to the Mott insulator, at which the 
correlation function for the bosonic superfluid Cb(t) = 
( c l h+r c b,h) switches from algebraic to exponential decay. 
We choose h = L/2. In Fig. [2] we plot Cb(V), the pola- 
ronic correlation function Cp(r) = (cb.h+rCj h+r c f>h c b u) 
and the polaron-pair correlation function Cpp(r) = 

<(Q£ h )toH, where °7 = cu+Aj+ivAj-, C ?(r) 
decays algebraically and oscillates with the Fermi wave- 
number kf = irvf/d as Cp(r) ~ sm(kfrd)\r\~ ap at long 
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FIG. 3: (color online) The exponents qb (blue squares), ap 
(green diamonds), app (red circles), and osdw (black trian- 
gles). While t/Ubb is varied at Ubf /Ubb = 1.4 in (a), Ubf /Ubb is 
varied at t/Ubb — 0.05 in (b). The red solid and green dashed 
lines represent app = 2/K p and ap = K p /2 + 1/(2K P ), 
which are derived from TLL theory for the ID polaron gas 
of Eq. © [1^]. The Luttinger parameter is given by K p = 
71-/(271- - 2cos _1 (l - u)) at u p = 0.5 [3C|. 




FIG. 4: (color online) The static structure factor for the spin 
density 5sd(A:) for Ubf /Ubb = 0.4 (blue solid) and 1.4 (red 
dashed). We set Ut = Vj = 0.5 and t/Ubb — 0.05 at which the 
system is deep in the Mott insulating state. 



distances. Cpp(r) = ((Of^ h )^O^ p ) decays algebraically 
as Cpp(r) - |r|- Qpp . As shown in Fig.^a) (Fig.^b)), 
C^(r) (Cp(r)) decays more slowly in the two-component 
TLL (Mott insulating) regime. By fitting the correlation 
functions with f(r) = 7|?'| _a , where a and 7 arc fitting 
parameters, we extract the exponents a and plot them 
as functions of t in Fig. Eta). For small i, the exponent 
for the bosonic correlation as exceeds that of the pola- 
ronic correlation ap at what we define as the transition 
point. By this, we locate the Mott transition as shown 
by the purple thick-solid lines in Figs.QJb) and (c). The 
numerical phase boundary, for L = 80, is closer to the 
strong-coupling regime determined by the RG analyses 
than the actual phase boundary, which is expected to 
emerge for larger systems. Note that t he p hase bound- 
ary agrees well with the results in Ref. [29( , depicted by 
the orange dashed-dotted line in Fig. [TJb). 

As discussed above, there are three phases in the Mott 
insulating regime, namely PS, PP-CFSF, and SDW, 
see [4o| . PS is signaled by an emerging peak in the 
bosonic structure factor at low wave-number (29| . The 
PP-CFSF phase is characterized by the correlation func- 
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FIG. 5: (color online) Ground state properties in the presence 
of a parabolic trapping potential, where Nt = Nf = 20, Q/t = 
0.008, and Ubf/Ubb — 1.4. (a) Density profiles pj for several 
values of t/Ubb are shown, (b) The correlation functions Cb(t) 
and Cpp(r) are plotted on a log-log scale, where t/Ubb = 0.05. 



tion Cpp(r). In Fig. [3ta), we see that ctpp < 2 in the 
Mott insulating regime at u = 1.4, meaning that the sys- 
tem has PP-CFSF QLRO. To identify SDW QLRO, we 
calculate the static structure factor for the spin density 



-ikd(j-l) 



(7) 



3,1 



where Arij = n&j — n f.j- Since (And +r Ani) ~ 
cos(2fc/re?)|r|~ QSDW at long distance, the structure factor 
behaves as Ssn(k) ~ ||fc| - 2fc / | QSDW - 1 near |fc| = 2kf. 
Consequently, Ssu(k) has cusps at k = ±2fc/ when the 
SDW QLRO is present, i.e. a SD w < 2. In Fig. gj we 
show Ssu(k) at v\, = Vj = 0.5 and indeed find a cusp for 
u = 0.4. By numerically fitting Ssd(&) near |fc| = 2kf, 
we extract the exponent asDW- I r Fig. [StL), agrjw is 
plotted as a function of u for i — 0.05, compared to app. 
Clearly the system transitions from SDW to PP-CFSF 
when u increases. As seen in Figs. QJb) an d ( c ), the 
PP-CFSF phase is present approximately for 1 < u < 2 
and a wide range of Ais, consistent with the analytical 
results above. Note that in previous studies the polaron- 
pair correlation function was not considered, and thus 
the PP-CFSF phase not found. 

Given that u = 1.32 in the 174 Yb- 173 Yb mixture of 
Ref. 0, PP-CFSF is expected to be realized in the dual 
Mott insulator regime when the temperature is suffi- 
ciently low. To demonstrate the stability of the phase 
in a trap, we calculate the ground states including a 



parabolic trap, ej = 



', where il and j c are the 



curvature and center of the trap. We set Q/t = 0.008 and 
Nb = Nf = 20, where N a is the number of particle a. In 
Fig. E^a), we show the density profiles pj = (rib,j + nfj) 
for several t while fixing u = 1.4. At t = 0.35, the system 
is compressible and no Mott plateau is present. When i 
decreases, a plateau at v = 1 is formed around i ~ 0.18. 
We calculate Cs(r) and Cpp(r) at t = 0.05 as shown 
in Fig. EJb). Inside the plateau, Cb(»") decays exponen- 
tially while Cpp(r) decays algebraically with the expo- 
nent a PP = 1.60 ± 0.05 < 2. Thus, there is PP-CFSF 



QLRO. 

Finally, we estimate the temperature regime required 
for PP-CFSF to emerge. At finite temperatures (T > 0) 
in ID, the polaron-pair correlation function Cpp(r) ex- 
hibits algebraic decay at r <C £ while decaying ex- 
ponentially for r ^> £, where the crossover length is 
£ = Vp/{ixk^,T) and the sound velocity of the polarons 
is v p = 2Jdsin(7Ti/ p )(l + V/(wJ)) for small V/J [U[39j]. 
Hence, when £ > d, PP-CFSF QLRO is considered to 
be sufficiently developed so that signatures of the PP- 
CFSF state can be observed. Since the condition £ ;§> d 
corresponds to J > k^T, we estimate J, as an example, 
for t/Ubb = 0.18, which is a modest value for a com- 
plete Mott plateau to be formed as shown in Fig. [5ja) . 
In the experiments of 174 Yb- 173 Yb mixtures, the lattice 
spacing is d = 266nm and the s-wave scattering lengths 
are abb = 5.55nm and a,bf = 7.34nm. Moreover, we as- 
sume the lattice depth in the transverse direction to be 
V± = 50-Er, where En is the recoil energy. Using these 
parameters and the Wannier function obtained by numer- 
ically solving the Schrbdingcr equation with a sinusoidal 
potential, we estimate that t/Ubb = 0.18 is reached when 
the lattice depth in the axial direction is Vb — 2.33£ , r. 
At this lattice depth, J ~ fee x 7.1nK. Since the lowest 
temperature realized in Ref. Q is ~ 5nK, the condition 
k^T < J is already possible with current experimental 
techniques. Further experimental advancement is neces- 
sary to achieve k-gT <C J ■ 

In summary, we have used the time-evolving block dec- 
imation method and Tomonaga-Luttinger liquid theory 
to reveal the quantum phases of polarons inside the dual 
Mott insulator of one-dimensional Bose-Fermi mixtures 
at unit total filling. Interestingly, we found a large phase 
diagram regime with a counterfiow superfluid phase of 
polaron pairs (PP-CFSF), in contrast to previous stud- 
ies. We have shown that this state is expected to be 
formed in the dual Mott insulator regime of 174 Yb- 173 Yb 
mixtures in optical lattices at sufficiently low tempera- 
tures. 
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